home *** CD-ROM | disk | FTP | other *** search
/ QRZ! Ham Radio 8 / QRZ Ham Radio Callsign Database - Volume 8.iso / mac / files / ant_nec / nec81tar.z / nec81tar / eksc.f < prev    next >
Text File  |  1991-05-13  |  14KB  |  579 lines
  1. C $TITLE: 'EKSC'
  2. C $NOFLOATCALLS
  3. C
  4.       SUBROUTINE EKSC (S,Z,RH,XK,IJ,EZS,ERS,EZC,ERC,EZK,ERK)
  5. C     COMPUTE E FIELD OF SINE, COSINE, AND CONSTANT CURRENT FILAMENTS BY
  6. C     THIN WIRE APPROXIMATION.
  7.       REAL*8 CONX(2),XK,RKB2,RHK,SHK,SS,CS
  8.       COMPLEX*16 CON,GZ1,GZ2,GP1,GP2,GZP1,GZP2,EZS,ERS,EZC,ERC,EZK,ERK
  9.       COMMON /TMI/ ZPK,RKB2,IJX
  10.       EQUIVALENCE (CONX,CON)
  11.       DATA CONX/0.,4.771341189D0/
  12. C**
  13.       IJX=IJ
  14.       ZPK=XK*Z
  15.       RHK=XK*RH
  16.       RKB2=RHK*RHK
  17.       SH=.5*S
  18.       SHK=XK*SH
  19.       SS=DSIN(SHK)
  20.       CS=DCOS(SHK)
  21.       Z2=SH-Z
  22.       Z1=-(SH+Z)
  23. C**
  24.       CALL GX (Z1,RH,XK,GZ1,GP1)
  25.       CALL GX (Z2,RH,XK,GZ2,GP2)
  26. C**
  27.       GZP1=GP1*Z1
  28.       GZP2=GP2*Z2
  29.       EZS=CON*((GZ2-GZ1)*CS*XK-(GZP2+GZP1)*SS)
  30.       EZC=-CON*((GZ2+GZ1)*SS*XK+(GZP2-GZP1)*CS)
  31.       ERK=CON*(GP2-GP1)*RH
  32.       CALL INTX (-SHK,SHK,RHK,IJ,CINT,SINT)
  33.       EZK=-CON*(GZP2-GZP1+XK*XK*CMPLX(CINT,-SINT))
  34.       GZP1=GZP1*Z1
  35.       GZP2=GZP2*Z2
  36.       IF (RH.LT.1.E-10) GO TO 1
  37.       ERS=-CON*((GZP2+GZP1+GZ2+GZ1)*SS-(Z2*GZ2-Z1*GZ1)*CS*XK)/RH
  38.       ERC=-CON*((GZP2-GZP1+GZ2-GZ1)*CS+(Z2*GZ2+Z1*GZ1)*SS*XK)/RH
  39. C**
  40.       RETURN
  41. 1     ERS=(0.,0.)
  42.       ERC=(0.,0.)
  43.       RETURN
  44.       END
  45. C
  46. C
  47. C
  48.       SUBROUTINE EKSCX(BX,S,Z,RHX,XK,IJ,INX1,INX2,EZS,ERS,EZC,ERC,EZK,
  49.      1ERK)
  50. C     COMPUTE E FIELD OF SINE, COSINE, AND CONSTANT CURRENT FILAMENTS BY
  51. C     EXTENDED THIN WIRE APPROXIMATION.
  52.       INTEGER*4 INX1,INX2
  53.       REAL*8 CONX(2),XK,RKB2,RHK,SHK,SS,CS
  54.       COMPLEX*16 CON,EZS,ERS,EZC,ERC,EZK,ERK,GZ1,GZ2,GRK1,GRK2
  55.       COMPLEX*16 GZP1,GZP2,GR1,GR2,GRP1,GRP2,GZZ1,GZZ2
  56.       COMMON /TMI/ ZPK,RKB2,IJX
  57.       EQUIVALENCE (CONX,CON)
  58.       DATA CONX/0.,4.771341189D0/
  59. C**
  60.       IF (RHX.LT.BX) GO TO 1
  61.       RH=RHX
  62.       B=BX
  63.       IRA=0
  64.       GO TO 2
  65. 1     RH=BX
  66.       B=RHX
  67.       IRA=1
  68. 2     SH=.5*S
  69.       IJX=IJ
  70.       ZPK=XK*Z
  71.       RHK=XK*RH
  72.       RKB2=RHK*RHK
  73.       SHK=XK*SH
  74.       SS=DSIN(1.D0*SHK)
  75.       CS=DCOS(1.D0*SHK)
  76.       Z2=SH-Z
  77.       Z1=-(SH+Z)
  78.       A2=B*B
  79.       IF (INX1.EQ.2) GO TO 3
  80.       CALL GXX (Z1,RH,B,A2,XK,IRA,GZ1,GZP1,GR1,GRP1,GRK1,GZZ1)
  81.       GO TO 4
  82. 3     CALL GX (Z1,RHX,XK,GZ1,GRK1)
  83.       GZP1=GRK1*Z1
  84.       GR1=GZ1/RHX
  85.       GRP1=GZP1/RHX
  86.       GRK1=GRK1*RHX
  87.       GZZ1=(0.,0.)
  88. 4     IF (INX2.EQ.2) GO TO 5
  89.       CALL GXX (Z2,RH,B,A2,XK,IRA,GZ2,GZP2,GR2,GRP2,GRK2,GZZ2)
  90.       GO TO 6
  91. 5     CALL GX (Z2,RHX,XK,GZ2,GRK2)
  92.       GZP2=GRK2*Z2
  93.       GR2=GZ2/RHX
  94.       GRP2=GZP2/RHX
  95.       GRK2=GRK2*RHX
  96.       GZZ2=(0.,0.)
  97. 6     EZS=CON*((GZ2-GZ1)*CS*XK-(GZP2+GZP1)*SS)
  98.       EZC=-CON*((GZ2+GZ1)*SS*XK+(GZP2-GZP1)*CS)
  99.       ERS=-CON*((Z2*GRP2+Z1*GRP1+GR2+GR1)*SS-(Z2*GR2-Z1*GR1)*CS*XK)
  100.       ERC=-CON*((Z2*GRP2-Z1*GRP1+GR2-GR1)*CS+(Z2*GR2+Z1*GR1)*SS*XK)
  101.       ERK=CON*(GRK2-GRK1)
  102.       CALL INTX (-SHK,SHK,RHK,IJ,CINT,SINT)
  103.       BK=B*XK
  104.       BK2=BK*BK*.25
  105.       EZK=-CON*(GZP2-GZP1+XK*XK*(1.-BK2)*CMPLX(CINT,-SINT)-BK2*(GZZ2
  106.      1 -GZZ1))
  107.       RETURN
  108.       END
  109. C
  110. C
  111. C
  112.       SUBROUTINE ROM2 (A,B,SUM,DMIN)
  113. C
  114. C     FOR THE SOMMERFELD GROUND OPTION, ROM2 INTEGRATES OVER THE SOURCE
  115. C     SEGMENT TO OBTAIN THE TOTAL FIELD DUE TO GROUND.  THE METHOD OF
  116. C     VARIABLE INTERVAL WIDTH ROMBERG INTEGRATION IS USED.  THERE ARE 9
  117. C     FIELD COMPONENTS - THE X, Y, AND Z COMPONENTS DUE TO CONSTANT,
  118. C     SINE, AND COSINE CURRENT DISTRIBUTIONS.
  119. C
  120.       REAL*8 A,B,TMAG1,TMAG2,TR,TI,DMIN,Z,DZ,DZOT,S,EP,ZE,ZEND
  121.       COMPLEX*16 SUM,G1,G2,G3,G4,G5,T00,T01,T10,T02,T11,T20
  122.       INTEGER*4 NM,NS
  123.       DIMENSION SUM(9),G1(9),G2(9),G3(9),G4(9),G5(9),T01(9),T10(9),
  124.      1T20(9)
  125.       DATA NM,NTS,NX,N/65536,4,1,9/,RX/1.E-4/
  126. C**
  127.       Z=A
  128.       ZE=B
  129.       S=B-A
  130.       IF (S.GE.0.) GO TO 1
  131.       WRITE(*,18)
  132.       STOP
  133. 1     EP=S/(1.D4*NM)
  134.       ZEND=ZE-EP
  135.       DO 2 I=1,N
  136. 2     SUM(I)=(0.,0.)
  137.       NS=NX
  138.       NT=0
  139.       CALL SFLDS (Z,G1)
  140. 3     DZ=S/NS
  141.       IF (Z+DZ.LE.ZE) GO TO 4
  142.       DZ=ZE-Z
  143.       IF (DZ.LE.EP) GO TO 17
  144. 4     DZOT=DZ*.5
  145.       CALL SFLDS (Z+DZOT,G3)
  146.       CALL SFLDS (Z+DZ,G5)
  147. 5     TMAG1=0.
  148.       TMAG2=0.
  149. C
  150. C     EVALUATE 3 POINT ROMBERG RESULT AND TEST CONVERGENCE.
  151. C
  152.       DO 6 I=1,N
  153.       T00=(G1(I)+G5(I))*DZOT
  154.       T01(I)=(T00+DZ*G3(I))*.5
  155.       T10(I)=(4.*T01(I)-T00)/3.
  156.       IF (I.GT.3) GO TO 6
  157.       TR=REAL(T01(I))
  158. C      TI=AIMAG(T01(I))
  159.       TI=IMAG(T01(I))
  160.       TMAG1=TMAG1+TR*TR+TI*TI
  161.       TR=REAL(T10(I))
  162. C      TI=AIMAG(T10(I))
  163.       TI=IMAG(T10(I))
  164.       TMAG2=TMAG2+TR*TR+TI*TI
  165. 6     CONTINUE
  166.       TMAG1=DSQRT(TMAG1)
  167.       TMAG2=DSQRT(TMAG2)
  168.       CALL TEST(TMAG1,TMAG2,TR,0.D0,0.D0,TI,DMIN)
  169.       IF(TR.GT.RX)GO TO 8
  170.       DO 7 I=1,N
  171. 7     SUM(I)=SUM(I)+T10(I)
  172.       NT=NT+2
  173.       GO TO 12
  174. 8     CALL SFLDS (Z+DZ*.25,G2)
  175.       CALL SFLDS (Z+DZ*.75,G4)
  176.       TMAG1=0.
  177.       TMAG2=0.
  178. C
  179. C     EVALUATE 5 POINT ROMBERG RESULT AND TEST CONVERGENCE.
  180. C
  181.       DO 9 I=1,N
  182.       T02=(T01(I)+DZOT*(G2(I)+G4(I)))*.5
  183.       T11=(4.*T02-T01(I))/3.
  184.       T20(I)=(16.*T11-T10(I))/15.
  185.       IF (I.GT.3) GO TO 9
  186.       TR=REAL(T11)
  187. C      TI=AIMAG(T11)
  188.       TI=IMAG(T11)
  189.       TMAG1=TMAG1+TR*TR+TI*TI
  190.       TR=REAL(T20(I))
  191. C      TI=AIMAG(T20(I))
  192.       TI=IMAG(T20(I))
  193.       TMAG2=TMAG2+TR*TR+TI*TI
  194. 9     CONTINUE
  195.       TMAG1=DSQRT(TMAG1)
  196.       TMAG2=DSQRT(TMAG2)
  197.       CALL TEST(TMAG1,TMAG2,TR,0.D0,0.D0,TI,DMIN)
  198.       IF(TR.GT.RX)GO TO 14
  199. 10    DO 11 I=1,N
  200. 11    SUM(I)=SUM(I)+T20(I)
  201.       NT=NT+1
  202. 12    Z=Z+DZ
  203.       IF (Z.GT.ZEND) GO TO 17
  204.       DO 13 I=1,N
  205. 13    G1(I)=G5(I)
  206.       IF (NT.LT.NTS.OR.NS.LE.NX) GO TO 3
  207.       NS=NS/2
  208.       NT=1
  209.       GO TO 3
  210. 14    NT=0
  211.       IF (NS.LT.NM) GO TO 15
  212.       WRITE(*,19)  Z
  213.       GO TO 10
  214. 15    NS=NS*2
  215.       DZ=S/NS
  216.       DZOT=DZ*.5
  217.       DO 16 I=1,N
  218.       G5(I)=G3(I)
  219. 16    G3(I)=G2(I)
  220.       GO TO 5
  221. 17    CONTINUE
  222. C**
  223.       RETURN
  224. C
  225. 18    FORMAT (' ERROR - B LESS THAN A IN ROM2')
  226. 19    FORMAT (33H ROM2 -- STEP SIZE LIMITED AT Z =,1P,E12.5)
  227.       END
  228. C
  229. C
  230. C
  231.       SUBROUTINE GX (ZZ,RH,XK,GZ,GZP)
  232. C     SEGMENT END CONTRIBUTIONS FOR THIN WIRE APPROX.
  233.       REAL*8 XK,R,R2,RK,CRK,SRK
  234.       COMPLEX*16 GZ,GZP
  235. C**
  236.       R2=ZZ*ZZ+RH*RH
  237.       R=DSQRT(R2)
  238.       RK=XK*R
  239.       CRK= DCOS(RK)
  240.       SRK=-DSIN(RK)
  241.       GZ=DCMPLX(CRK,SRK)/R
  242.       GZP=-DCMPLX(1.,RK)*GZ/R2
  243.       END
  244. C
  245. C
  246. C
  247.       SUBROUTINE GXX (ZZ,RH,A,A2,XK,IRA,G1,G1P,G2,G2P,G3,GZP)
  248. C     SEGMENT END CONTRIBUTIONS FOR EXT. THIN WIRE APPROX.
  249.       REAL*8 XK,R,R2,R4,RK,RK2,RH2
  250.       COMPLEX*16 GZ,C1,C2,C3,G1,G1P,G2,G2P,G3,GZP
  251.       R2=ZZ*ZZ+RH*RH
  252.       R=DSQRT(R2)
  253.       R4=R2*R2
  254.       RK=XK*R
  255.       RK2=RK*RK
  256.       RH2=RH*RH
  257.       T1=.25*A2*RH2/R4
  258.       T2=.5*A2/R2
  259.       C1=CMPLX(1.,RK)
  260.       C2=3.*C1-RK2
  261.       C3=CMPLX(6.,RK)*RK2-15.*C1
  262.       GZ=CMPLX(DCOS(RK),-DSIN(RK))/R
  263.       G2=GZ*(1.+T1*C2)
  264.       G1=G2-T2*C1*GZ
  265.       GZ=GZ/R2
  266.       G2P=GZ*(T1*C3-C1)
  267.       GZP=T2*C2*GZ
  268.       G3=G2P+GZP
  269.       G1P=G3*ZZ
  270.       IF (IRA.EQ.1) GO TO 2
  271.       G3=(G3+GZP)*RH
  272.       GZP=-ZZ*C1*GZ
  273.       IF (RH.GT.1.E-10) GO TO 1
  274.       G2=0.
  275.       G2P=0.
  276.       RETURN
  277. 1     G2=G2/RH
  278.       G2P=G2P*ZZ/RH
  279.       RETURN
  280. 2     T2=.5*A
  281.       G2=-T2*C1*GZ
  282.       G2P=T2*GZ*C2/R2
  283.       G3=RH2*G2P-A*GZ*C1
  284.       G2P=G2P*ZZ
  285.       GZP=-ZZ*C1*GZ
  286.       RETURN
  287.       END
  288. C
  289. C
  290. C
  291.       SUBROUTINE INTX (EL1,EL2,B,IJ,SGR,SGI)
  292. C
  293. C     INTX PERFORMS NUMERICAL INTEGRATION OF EXP(JKR)/R BY THE METHOD OF
  294. C     VARIABLE INTERVAL WIDTH ROMBERG INTEGRATION.  THE INTEGRAND VALUE
  295. C     IS SUPPLIED BY SUBROUTINE GF.
  296. C
  297.       INTEGER*4 NM,NS
  298.       REAL*8 T01R,T10R,TE1R,T01I,T10I,TE1I,T11R,T20R,TE2R,TE2I,T11I,
  299.      1 T20I,TE21,T00R,T00I,G1R,G1I,G2R,G2I,G3R,G3I,G4R,G4I,G5R,G5I,
  300.      2 RKB2,B,S,EL1,EL2,Z,ZE,ZP,DZ,DZOT
  301.       COMMON /TMI/ ZPK,RKB2,IJX
  302.       DATA NX,NM,NTS,RX/1,65536,4,1.E-4/
  303. C**
  304.       Z=EL1
  305.       ZE=EL2
  306.       IF (IJ.EQ.0) ZE=0.
  307.       S=ZE-Z
  308.       FNM=NM
  309.       EP=S/(10.*FNM)
  310.       ZEND=ZE-EP
  311.       SGR=0.
  312.       SGI=0.
  313.       NS=NX
  314.       NT=0
  315. C**
  316.       CALL GF (Z,G1R,G1I)
  317. C**
  318. 1     FNS=NS
  319.       DZ=S/FNS
  320.       ZP=Z+DZ
  321.       IF (ZP-ZE) 3,3,2
  322. 2     DZ=ZE-Z
  323.       IF (ABS(DZ)-EP) 17,17,3
  324. 3     DZOT=DZ*.5
  325.       ZP=Z+DZOT
  326.       CALL GF (ZP,G3R,G3I)
  327.       ZP=Z+DZ
  328.       CALL GF (ZP,G5R,G5I)
  329. 4     T00R=(G1R+G5R)*DZOT
  330.       T00I=(G1I+G5I)*DZOT
  331.       T01R=(T00R+DZ*G3R)*0.5
  332.       T01I=(T00I+DZ*G3I)*0.5
  333.       T10R=(4.0*T01R-T00R)/3.0
  334.       T10I=(4.0*T01I-T00I)/3.0
  335. C
  336. C     TEST CONVERGENCE OF 3 POINT ROMBERG RESULT.
  337. C
  338.       CALL TEST (T01R,T10R,TE1R,T01I,T10I,TE1I,0.D0)
  339.       IF (TE1I-RX) 5,5,6
  340. 5     IF (TE1R-RX) 8,8,6
  341. 6     ZP=Z+DZ*0.25
  342.       CALL GF (ZP,G2R,G2I)
  343.       ZP=Z+DZ*0.75
  344.       CALL GF (ZP,G4R,G4I)
  345.       T02R=(T01R+DZOT*(G2R+G4R))*0.5
  346.       T02I=(T01I+DZOT*(G2I+G4I))*0.5
  347.       T11R=(4.0*T02R-T01R)/3.0
  348.       T11I=(4.0*T02I-T01I)/3.0
  349.       T20R=(16.0*T11R-T10R)/15.0
  350.       T20I=(16.0*T11I-T10I)/15.0
  351. C
  352. C     TEST CONVERGENCE OF 5 POINT ROMBERG RESULT.
  353. C
  354.       CALL TEST (T11R,T20R,TE2R,T11I,T20I,TE2I,0.D0)
  355.       IF (TE2I-RX) 7,7,14
  356. 7     IF (TE2R-RX) 9,9,14
  357. 8     SGR=SGR+T10R
  358.       SGI=SGI+T10I
  359.       NT=NT+2
  360.       GO TO 10
  361. 9     SGR=SGR+T20R
  362.       SGI=SGI+T20I
  363.       NT=NT+1
  364. 10    Z=Z+DZ
  365.       IF (Z-ZEND) 11,17,17
  366. 11    G1R=G5R
  367.       G1I=G5I
  368.       IF (NT-NTS) 1,12,12
  369. 12    IF (NS-NX) 1,1,13
  370. C
  371. C     DOUBLE STEP SIZE
  372. C
  373. 13    NS=NS/2
  374.       NT=1
  375.       GO TO 1
  376. 14    NT=0
  377.       IF (NS-NM) 16,15,15
  378. 15    WRITE(*,20)  Z
  379.       GO TO 9
  380. C
  381. C     HALVE STEP SIZE
  382. C
  383. 16    NS=NS*2
  384.       FNS=NS
  385.       DZ=S/FNS
  386.       DZOT=DZ*0.5
  387.       G5R=G3R
  388.       G5I=G3I
  389.       G3R=G2R
  390.       G3I=G2I
  391.       GO TO 4
  392. 17    CONTINUE
  393.       IF (IJ) 19,18,19
  394. C
  395. C     ADD CONTRIBUTION OF NEAR SINGULARITY FOR DIAGONAL TERM
  396. C
  397. 18    SGR=2.*(SGR+DLOG((DSQRT(B*B+S*S)+S)/B))
  398.       SGI=2.*SGI
  399. 19    CONTINUE
  400.       RETURN
  401. C
  402. 20    FORMAT (' STEP SIZE LIMITED AT Z=',F10.5)
  403.       END
  404. C
  405. C
  406. C
  407.       SUBROUTINE GF (ZK,CO,SI)
  408. C
  409. C     GF COMPUTES THE INTEGRAND EXP(JKR)/(KR) FOR NUMERICAL INTEGRATION.
  410. C
  411.       REAL*8 ZK,CO,SI,RK,RKS,RKB2
  412.       COMMON /TMI/ ZPK,RKB2,IJ
  413. C**
  414.       ZDK=ZK-ZPK
  415.       RK=DSQRT(RKB2+ZDK*ZDK)
  416.       SI=DSIN(RK)/RK
  417.       IF (IJ) 1,2,1
  418. 1     CO=DCOS(RK)/RK
  419.       RETURN
  420. 2     IF (RK.LT..2D0) GO TO 3
  421.       CO=(DCOS(RK)-1.)/RK
  422.       RETURN
  423. 3     RKS=RK*RK
  424.       CO=((-1.38888889D-3*RKS+4.16666667D-2)*RKS-.5)*RK
  425.       RETURN
  426.       END
  427. C
  428. C
  429. C
  430.       SUBROUTINE TEST (F1R,F2R,TR,F1I,F2I,TI,DMIN)
  431. C
  432. C     TEST FOR CONVERGENCE IN NUMERICAL INTEGRATION
  433. C
  434.       REAL*8 F1R,F2R,TR,F1I,F2I,TI,DMIN,DEN
  435.       DEN=ABS(F2R)
  436.       TR=ABS(F2I)
  437.       IF (DEN.LT.TR) DEN=TR
  438.       IF (DEN.LT.DMIN) DEN=DMIN
  439.       IF (DEN.LT.1.D-37) GO TO 1
  440.       TR=ABS((F1R-F2R)/DEN)
  441.       TI=ABS((F1I-F2I)/DEN)
  442.       RETURN
  443. 1     TR=0.
  444.       TI=0.
  445.       RETURN
  446.       END
  447. C
  448. C
  449. C
  450.       SUBROUTINE SFLDS (T,E)
  451. C
  452. C     SFLDX RETURNS THE FIELD DUE TO GROUND FOR A CURRENT ELEMENT ON
  453. C     THE SOURCE SEGMENT AT T RELATIVE TO THE SEGMENT CENTER.
  454. C
  455.       REAL*8 T,PI,TP,THET,POT,RK,R1,R2,ZMH,ZPH,SFAC
  456.       COMPLEX*16 E,U,U2,XX1,XX2,ERV,EZV,ERH,EZH,EPH
  457.       COMPLEX*16 T1,EXK,EYK,EZK,EXS,EYS,EZS,EXC,EYC,EZC,
  458.      1 ZRATI,ZRATI2,FRATI,ER,ET,HRV,HZV,HRH
  459.       INTEGER*4 IND1,IND2
  460.       COMMON/DATAJ/S,B,XJ,YJ,ZJ,CABJ,SABJ,SALPJ,EXK,EYK,EZK,EXS,EYS,
  461.      1 EZS,EXC,EYC,EZC,RKH,IND1,IND2,IPGND,IEXK
  462.       COMMON /INCOM/ XO,YO,ZO,SN,XSN,YSN,ISNOR
  463.       COMMON /GWAV/ U,U2,XX1,XX2,R1,R2,ZMH,ZPH
  464.       COMMON /GND/ZRATI,ZRATI2,FRATI,CL,CH,SCRWL,SCRWR,NRADL,KSYMP,
  465.      1 IFAR,IPERF,T1,T2
  466.       DIMENSION E(9)
  467.       DATA PI/3.141592654D0/,TP/6.283185308D0/,POT/1.570796327D0/
  468. C**
  469. C     E      WRITE(*,*) '   SFLDS: START'
  470. C**
  471.       XT=XJ+T*CABJ
  472.       YT=YJ+T*SABJ
  473.       ZT=ZJ+T*SALPJ
  474.       RHX=XO-XT
  475.       RHY=YO-YT
  476.       RHS=RHX*RHX+RHY*RHY
  477. C      RHO=DSQRT(RHS)
  478.       RHO=SQRT(RHS)
  479.       IF (RHO.GT.0.) GO TO 1
  480.       RHX=1.
  481.       RHY=0.
  482.       PHX=0.
  483.       PHY=1.
  484.       GO TO 2
  485. 1     RHX=RHX/RHO
  486.       RHY=RHY/RHO
  487.       PHX=-RHY
  488.       PHY=RHX
  489. 2     CPH=RHX*XSN+RHY*YSN
  490.       SPH=RHY*XSN-RHX*YSN
  491.       IF (ABS(CPH).LT.1.E-10) CPH=0.
  492.       IF (ABS(SPH).LT.1.E-10) SPH=0.
  493.       ZPH=ZO+ZT
  494.       ZPHS=ZPH*ZPH
  495.       R2S=RHS+ZPHS
  496. C      R2=DSQRT(R2S)
  497.       R2=SQRT(R2S)
  498.       RK=R2*TP
  499.       XX2=CMPLX(DCOS(RK),-DSIN(RK))
  500.       IF (ISNOR.EQ.1) GO TO 3
  501. C
  502. C     USE NORTON APPROXIMATION FOR FIELD DUE TO GROUND.  CURRENT IS
  503. C     LUMPED AT SEGMENT CENTER WITH CURRENT MOMENT FOR CONSTANT, SINE,
  504. C     OR COSINE DISTRIBUTION.
  505. C
  506.       ZMH=1.
  507.       R1=1.
  508.       XX1=0.
  509. C**
  510.       CALL GWAVE (ERV,EZV,ERH,EZH,EPH)
  511. C**
  512.       ET=-(0.,4.77134)*FRATI*XX2/(R2S*R2)
  513.       ER=2.*ET*CMPLX(1.,RK)
  514.       ET=ET*CMPLX(1.-RK*RK,RK)
  515.       HRV=(ER+ET)*RHO*ZPH/R2S
  516.       HZV=(ZPHS*ER-RHS*ET)/R2S
  517.       HRH=(RHS*ER-ZPHS*ET)/R2S
  518.       ERV=ERV-HRV
  519.       EZV=EZV-HZV
  520.       ERH=ERH+HRH
  521.       EZH=EZH+HRV
  522.       EPH=EPH+ET
  523.       ERV=ERV*SALPJ
  524.       EZV=EZV*SALPJ
  525.       ERH=ERH*SN*CPH
  526.       EZH=EZH*SN*CPH
  527.       EPH=EPH*SN*SPH
  528.       ERH=ERV+ERH
  529.       E(1)=(ERH*RHX+EPH*PHX)*S
  530.       E(2)=(ERH*RHY+EPH*PHY)*S
  531.       E(3)=(EZV+EZH)*S
  532.       E(4)=0.
  533.       E(5)=0.
  534.       E(6)=0.
  535.       SFAC=PI*S
  536.       SFAC=DSIN(SFAC)/SFAC
  537.       E(7)=E(1)*SFAC
  538.       E(8)=E(2)*SFAC
  539.       E(9)=E(3)*SFAC
  540. C**
  541. C     E      WRITE(*,*) '   SFLDS: RETURN LINE 91'
  542. C**
  543.       RETURN
  544. C
  545. C     INTERPOLATE IN SOMMERFELD FIELD TABLES
  546. C
  547. 3     IF (RHO.LT.1.E-12) GO TO 4
  548.       THET=DATAN(ZPH/RHO)
  549.       GO TO 5
  550. 4     THET=POT
  551. 5     CALL INTRP (R2,THET,ERV,EZV,ERH,EPH)
  552. C     COMBINE VERTICAL AND HORIZONTAL COMPONENTS AND CONVERT TO X,Y,Z
  553. C     COMPONENTS.  MULTIPLY BY EXP(-JKR)/R.
  554.       XX2=XX2/R2
  555.       SFAC=SN*CPH
  556.       ERH=XX2*(SALPJ*ERV+SFAC*ERH)
  557.       EZH=XX2*(SALPJ*EZV-SFAC*ERV)
  558.       EPH=SN*SPH*XX2*EPH
  559. C     X,Y,Z FIELDS FOR CONSTANT CURRENT
  560.       E(1)=ERH*RHX+EPH*PHX
  561.       E(2)=ERH*RHY+EPH*PHY
  562.       E(3)=EZH
  563.       RK=TP*T
  564. C     X,Y,Z FIELDS FOR SINE CURRENT
  565.       SFAC=DSIN(RK)
  566.       E(4)=E(1)*SFAC
  567.       E(5)=E(2)*SFAC
  568.       E(6)=E(3)*SFAC
  569. C     X,Y,Z FIELDS FOR COSINE CURRENT
  570.       SFAC=DCOS(RK)
  571.       E(7)=E(1)*SFAC
  572.       E(8)=E(2)*SFAC
  573.       E(9)=E(3)*SFAC
  574. C**
  575. C     E      WRITE(*,*) '   SFLDS: RETURN LINE 125'
  576. C**
  577.       RETURN
  578.       END
  579.